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A simple nonlinear transmission-line model of the cochlea with longitudinal coupling is introduced 
that can reproduce Basilar membrane response and neural tuning in the chinchilla. It is found that 
the middle ear has little effect on cochlear resonances, and hence conclude that the theory of coherent 
reflections is not applicable to the model. The model also provides an explanation of the emergence 
of spontaneous otoacoustic emissions (SOAEs). It is argued that SOAEs arise from Hopf bifurcations 
of the transmission-line model and not from localized instabilities. The paper shows that emissions 
can become chaotic, intermittent and fragile to perturbations. 
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Significance: The cochlea is a remarkable device that 
out-performs any human-made system; it is sensitive to 
sounds over a million-fold intensity and a ten-octave fre¬ 
quency range, and can distinguish signals separated by 
microseconds at frequencies only 0.2% apart. Here we 
study the mechanisms that make this work. We present 
a nonlinear mathematical model that combines the key 
physiological processes, including both longitudinal cou¬ 
pling and hair cell motility, which produces response pat¬ 
terns that agree with experiments in different animals. A 
dynamical systems analysis of the model allows us chal¬ 
lenge existing theories on the source of spontaneous otoa¬ 
coustic emissions, suggesting that the entire organ, rather 
than localized instabilities, are key. 

The mammalian hearing organ is a sensitive sensory 
device that operates at the extremes of physical limits. 
It is capable of resolving sound pressure levels just above 
atmospheric thermal noise and of discriminating frequen¬ 
cies 0.2% apart m- In order to achieve these features 
the inner ear employs an active feedback mechanism [S], 
Like any feedback loop, it is possible for the one in the in¬ 
ner ear to become unstable. In this paper we show how 
such an instability can lead to self-excited oscillations 
that are emitted from the ear as sound, which we pro¬ 
pose as a mechanism for the generation of spontaneous 
otoacoustic emissions (SOAEs) [13]. We derive this from 
a new mathematical model of the mammalian ear, which 
includes both active somatic motility and longitudinal 
coupling in the cochlea, linked to the atmosphere via the 
middle ear. 

One widespread explanation of emission generation is 
due to Shera |22| . Shera and Zweig [33 and Talmadge 
et al. [28] , They argue that spontaneous emissions are 
wave instabilities that arise because of coherent reflec¬ 
tion at the stapes and at the characteristic frequency 
(CF) position of the cochlea; a mechanism that gives rise 
to standing waves that are stabilized by the cochlear non- 
linearities. However, this theory is based on linear prop- 
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erties of the cochlea and does not explain how standing 
waves are stabilized. 

Here, we apply the theory of dynamical systems to 
revisit the question of SOAEs, which leads to an alterna¬ 
tive explanation for emission generation. We show that 
stability loss in the model generically corresponds to a 
Hopf bifurcation El, which can lead to a sustained peri¬ 
odic emission. This is different from the so-called ‘Hopf 
ear’ [5], in that it is an emergent macro-level property 
of the cochlear model, with all its many interconnected 
elements, rather than being localized in any individual 
micro-scale component. Indeed, one advantage of our 
approach is that there is no need to assume, a priori , the 
mechanism of the stability loss; bifurcation theory deter¬ 
mines the onset of instability irrespective of its source. 
We further predict that these bifurcations are ubiqui¬ 
tous, provided there is fine-scale spatial variation of the 
cochlear parameters. This finding provides a consistent 
explanation for widespread observation of SOAEs; in re¬ 
alistic operating conditions a bifurcation point is never 
far away. As one goes beyond a Hopf bifurcation, the re¬ 
sulting small amplitude limit-cycle vibrations can grow 
in amplitude, deform, and undergo further instabilities 
and bifurcations. One possible physical manifestation of 
such complex motions is in the way emissions can ap¬ 
pear, disappear or change their character as the cochlea 
undergoes changes due to damage [32] or aging [T5], for 
example. 

While linear cochlear models can also predict instabil¬ 
ity, they cannot adequately explain what happens when 
multiple instabilities produce a multiple independent pe¬ 
riodic motions. In a linear system, multiple oscillations 
can coexist without any effect on each other. In nonlinear 
systems, such as the model we present here, we expect 
more complex phenomena to arise, such as intermittent 
or chaotic oscillations to arise, for which there is some 
experimental evidence. For example, Burns |3_ describes 
experimentally observed short-term property changes to 
SOAEs, including their seemingly random sudden ap¬ 
pearance and disappearance. Keefe et al. m used time 
series analysis methods that could distinguish chaotic 
spontaneous as well as stimulated emissions. It is though 
perhaps not surprising that there are not more reports of 
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intermittent or chaotic SOAEs. To eliminate environ¬ 
mental noise one must use temporal averaging. However, 
averaging can also mask short term spectral variations 
and chaos. 

The rest of this paper is organized as follows. First, we 
introduce our mathematical model. Then we show how 
the model reproduces chinchilla BM mechanical data and 
neural tuning. Applying spatial roughness to the feed¬ 
back coefficients of the model, we show that the cochlea 
develops resonances. As we increase the size of the rough¬ 
ness the resonances turn into instabilities and we observe 
spontaneous oscillations. We investigate these vibrations 
in detail and conclude that our model can produce com¬ 
plicated dynamics in line with experimental observations. 


I. NONLINEAR TRANSMISSION-LINE MODEL 
OF THE INNER EAR 


We start from a simple transmission-line model of the 
cochlea inspired by the work of Zweig |29j . Our model, 
instead of a delayed feedback, contains a distributed feed¬ 
forward mechanism to account for the organ of Corti 
(OOC) structure. This type of active mechanism was 
introduced by Geisler [B] and Steele and Lim |23J and suc¬ 
cessfully applied in many models. Moreover, parameters 
of such a model are easily interpreted in terms of the lon¬ 
gitudinal coupling within the OOC |31| . In [27] together 
with Epp, we showed how such spatial coupling can be 
directly compared with time-delayed feedback, but with 
better stability properties. 

As is common in modeling the mechanics of the 
cochlea, we assume that the fluid in the perilymph is in¬ 
compressible and inviscid, which yields a wave equation 
for the pressure difference p between the scala tympani 
and vestibuli, and the Basilar membrane (BM) displace¬ 
ment Sbm, 


*bm(£,t) = I 1 ) 

where dot denotes (partial) differentiation with respect to 
time t. We nondimensionalise the longitudinal distance 
£, so that the length of the cochlea becomes unity, 0 ^ 
£ ^ 1. The lumped parameter e is a function of the 
geometry of the cochlear chambers and the density of 
the perilymph fluid, and m is the mass surface density of 
the BM. 

Our model of the BM motion, and its amplification by 
outer hair cells (OHCs), is a simplified version of that de¬ 
veloped by O Maoleidigh and Jiilicher m- This model 
accounts for the active nonlinearity of hair bundles, cou¬ 
pled to hair cell elongation, for a detailed derivation see 
[25] . To that model we add longitudinal coupling, as de¬ 
scribed in E23- Specifically, the BM displacement Xb m 
and OHC charge q , at longitudinal distance £ from the 


stapes, are defined by the equations 
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The BM motion has natural frequency ojq and relative 
damping £, both of which we assume to depend on posi¬ 
tion £. The BM is forced both by the cochlear pressure 
difference p at position £, and also by the OHCs which 
exert a force due to the unique protein called prestin |T] 
embedded in its lateral wall. We assume that the effect of 
the OHCs is longitudinally distributed, due to the push- 
pull mechanism of the OOC EBEU, giving rise to the 
integral term on the right-hand side of (J2|. The longitu¬ 
dinal convolution kernel w(x) is assumed to have a strong 
positive feed-forward and a weak negative feed-backward 
component, 


w(x) = e 5g ' > sin(a; — S a ), 


and a characteristic feed-forward distance is given by 
L(£) = that is about twice as long at the apex 

than at the base. Note that the kernal w resembles a 
strongly damped wave and therefore might alternatively 
represent a secondary transmission line generated in the 
tunnel of Corti HI] or on the tectoral membrane [7]. 

The charge q inside the OHC is controlled by the me¬ 
chanically sensitive ion channels of the hair bundle [9]. 
Equation |3| models the capacitor of the OHC that in¬ 
tegrates ion currents flowing through the ion channels. 
The OHC has an active mechanism to regain its resting 
potential, which is modeled by the rate constant 7 . We 
assume that the open probability of the ion channels is 
described by a second-order Boltzmann function [18| 
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The piezoelectric property of the OHC is represented by 
the constant v. 

The boundary condition of the model at the apex is 
determined by the helicotrema, which we assume has no 
resistance to fluid motion, hence p(l,t) = 0. The bound¬ 
ary condition at the base is determined by the action of 
the stapes on the oval window. The stapes is in turn 
controlled by the dynamics of the middle ear which we 
model as a single-degree-of-freedom oscillator 

^ow + 2 £„ w .x' ow T u,\ m x ovr — m ow (p (0, /.) T Gp e ), (4) 


where p e is the sound pressure acting on the eardrum and 
the subscript l ow’ refers to the oval window. Hence the 
basal boundary condition can be determined from force 
balance to be 


dp 

9£ 


( 0 ,*) mk ow x ow (t). 
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A. Parameter fitting 

The model 0-0 is only complete once parameters 
have been specified; a full list is given in Table [TJ We fit¬ 
ted the model to the frequency response function (FRF) 
of the BM at 0 dB sound pressure level (SPL) of ani¬ 
mal N92 in the data of Rhode m- The sharpness of 
tuning was fitting simple exponential functions to £(£), 
L(£) and y(£) and by finding a suitable value for f q . The 
shape of the FRF was tuned by varying the parameters 
of w(h), that is /3, S g and S 0 . However, tuning the model 
at a single position is not sufficient, because we want the 
sharpness of tuning to be accurate at every BM position. 
Unfortunately this type of data from a single cochlea does 
not currently exist. Reasoning that BM tuning and neu¬ 
ral tuning are close, we instead use the neural tuning data 
of Recio-Spinoso et al. |2Dj to further fit the parameters 
£(£), L(£) and y(£) of our model. 

The nonlinear response of the model was tuned at the 
CF position of 6.6 kHz for the 6.6 kHz stimulus. Only 
the parameters a\, b±, a 2 and &2 of the open probability 
function Po were altered, while keeping the derivative at 
the equilibrium Pq( 0) constant, so that the 0 dB FRF 
persisted. The sensitivity of the model is then tuned by 
adjusting only two parameters: G, which scales the input 
signal, and A, which tunes the sensitivity of the hair 
bundles, and effectively determines the BM amplitude 
range. 

II. RESULTS & DISCUSSION 
A. BM response 

In order to find the BM response, we solve our model 
Q-Q numerically; see the Materials and Methods sec¬ 
tion for details. The results are presented in Fig.[l]in the 
form of thick colored lines. The superimposed thin black 
curves show the corresponding experimental data. 

Results at the CF place of 6.6 kHz, for which our model 
is tuned, are shown in Figure [TJa-d). Figures 0 ^) and 
(b) show input-output functions; BM displacement as 
a function of sound pressure level (SPL) at various fre¬ 
quencies, and as a function of SPL at various frequencies, 
respectively. The data at lower than CF frequencies are 
accurately matched by the results, at higher than CF 
frequencies the calculated response is more compressed 
at lower SPLs than the data. The maximal compres¬ 
sion rate, however, matches the data well. Figure [ljb) 
shows the frequency response functions at different SPLs. 
Again, the agreement with data is good. However, for 
higher amplitude stimuli the peak in the data widens al¬ 
most symmetrically, while our calculations show almost 
the same width; this is the same disagreement as in Fig¬ 
ure [lja) . One of the reasons why the peak of the tuning 
curve in our calculations does not broaden is because our 
convolution kernal w does not depend on SPL. One way 
to improve the model would be to allow the parameters of 


w to vary with p e . Possible physical causes of such vari¬ 
ation include tectoral membrane waves [ 7 . and tunnel of 
Corti flow 1111 - 

Sound compression can be measured by the rate of 
growth (RoG), defined as the slope of the log-log graph 
of BM displacement versus SPL. Figure [TJc) shows the 
RoG, as a function of frequency, for a range of different 
stimulus intensities. Our model shows good agreement 
with the data. As expected, compression is minimal for 
frequencies significantly lower than CF, then increases (so 
that the RoG can be as low as 0.1). Note that, due to 
interference, the RoG fluctuates significantly at certain 
frequencies especially above CF, which the model also 
reproduces. The phase variation of the BM response is 
illustrated in Figs. 0 d). Again we see good agreement 
between model predictions and experimental data. 

Even though our model has been fitted to data from a 
specific animal whose BM vibration was measured at the 
6.6 kHz CF position, we can compare its predictions to 
another dataset at a different frequency, measured in a 
different animal of the same species. Such data are shown 
in Figure[lje-h). The model is adjusted by changing only 
two parameters: the input sensitivity G and hair bundle 
sensitivity A. The numerical results are still reasonably 
close to the data and all the qualitative features are pre¬ 
served. The most obvious differences are that our model 
shows slightly sharper tuning than the animal, and that 
the phase variation along the frequency axis is somewhat 
greater in the experimental data (again possibly a result 
of the assumption that w is independent of stimulus). 


B. Stability of ideal and rough cochlea models 

Having shown the validity of our model, we now use 
it to investigate the stability of the cochlea. Dynamical 
systems theory states that if the spectrum of the system 
at an equilibrium state has only negative real parts, the 
system is stable and small perturbations to that equi¬ 
librium will decay. Since our cochlea model 0-0 is 
a partial differential equation with spatial feedback and 
feed-forward, its spectrum will typically contain a con¬ 
tinuous curve in the complex plane, spanning frequencies 
from the lowest to the highest audible frequencies. How¬ 
ever, our use of a numerical discretization scheme (which 
approximates the system by a large set of ordinary differ¬ 
ential equations) means that we can in practice compute 
the spectrum using a standard eigenvalue solver, result¬ 
ing in a discrete approximation to the spectrum, with a 
large number of discrete eigenvalues. 

Figure [2|a) shows computed spectra of the cochlea 
model, in the absence of stimulus. The imaginary axis 
is shown as a dashed line; any eigenvalue A above this 
line indicates instability of the cochlea. The green cir¬ 
cles in Figure [2ja) (mostly organized in the thick line in 
the middle of the graph) represent the spectrum of the 
cochlea described above; they show that Im(A) < 0 for 
all A, so the cochlea is stable to small disturbances. 
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TABLE I. Parameter values of the mathematical model 
Parameter Value Description 


£ 

0.0594 

m 

0.0055 g-cm" 2 

OJQ 

20.832e _4 ' 8354e - 0.1455 

c 

0.11e 14f 

u 

0.2221 

L 

0.005e°' 7? 

7 

12.5 x 5 -5 ms -1 

V 

0.72 nC 

Ihb 

0.390323 nA 

ai 

19.9873 

a 2 

16.3928 

bi 

0.0692021 

£>2 

0.0369987 

77-ow 

1.85 g • cm -2 

/now 

1100 cm -1 

Cow 

0.265 

P 

0.78 

S g 

0.74 

So 

-0.04 

Lo 

0.005 

hL 

0.7 


cochlear geometry constant 

mass density of the Basilar membrane 

undamped natural frequency of the Basilar membrane 

relative damping of the Basilar membrane 

gain of the outer hair cell 

characteristic feed-forward distance 

RC time constant of the outer hair cell 

piezoelectric constant of the outer hair cell 

maximum transduction current through the hair bundles 

parameter of Po 

parameter of Po 

parameter of Po 

parameter of Po 

surface density of the stapes 

coupling of oval window to perilymph 

relative damping of the middle ear 

convolution kernel Gaussian characteristic width 

convolution kernel Gaussian shift 

convolution kernel wave shift 

feed-forward distance at base 

rate of growth of the feed-forward distance 


This calculation is for an ideal cochlea; a real organ will 
in general have small geometric variations of its proper¬ 
ties along its length. In order to model such inhomo¬ 
geneities, we allow the feed-forward parameter f q to vary 
randomly about its notional value at each position of the 
cochlea. We assume this variation to be normally dis¬ 
tributed with zero mean and constant variance a 2 . Such 
randomness has been shown by Ku et al. m to induce 
reflections from the Basilar membrane back towards the 
stapes. In what follows we shall consider only one realiza¬ 
tion of the distribution, but we shall allow the variation 
to be scaled, in effect varying the standard deviation a. 
We shall refer to such a model the rough cochlea, with a 
a measure of the degree of roughness. 

The spectrum of one particular rough cochlea (with 
<7 = 0.06) is represented in Figure [2ja) by red x signs. 
We see that the continuous spectrum breaks up; the 
eigenvalues are scattered, and discrete eigenvalues jump 
out of the curve found for the smooth cochlea, some¬ 
times by a significant distance. The closer an eigenvalue 
gets to the imaginary axis, the more the cochlea becomes 
sensitive to a disturbance at the corresponding frequency, 
which lowers the hearing threshold at that frequency. If 
an eigenvalue crosses the imaginary axis, a Hopf bifur¬ 
cation occurs; as a result the cochlea becomes unstable 
and we would expect a spontaneous oscillation in the ab¬ 
sence of any stimulus. This explains why SOAEs occur at 
the frequencies where the hearing threshold has notches. 


We also see in Figure |2ja) that spectral points close to 
the imaginary axis are roughly periodically spaced in fre¬ 
quency. This agrees with data that spontaneous and 
stimulus frequency OAEs and the hearing threshold mi¬ 
crostructures are roughly periodic in frequency [14^. 


C. Coherent reflections 

The observation from Figure [2]( a) that spectral points 
close to the imaginary axis are regularly spaced in fre¬ 
quency might be explained by the theory of Zweig and 
Shera [32}, whose assumption is that resonances build 
up in the cochlea between the oval window and the CF 
position of the resonance. To test such a hypothesis we 
removed the middle ear from our model altogether, and 
made sure that there is no reflection from the oval win¬ 
dow by setting p(0,t) = 0. Calculating the stability of 
this model, with the same roughness as before, yields a 
remarkably similar spectrum (illustrated by the blue + 
signs in Figure |2ja)). This means that, while reflections 
from the oval window play a part in shaping resonances, 
they are not essential. 

To further clarify how reflections occur we also cal¬ 
culated the eigenvector that corresponds to the largest 
Re(A)/|Im(A)| value; shown in Figure §b ,c), with and 
without the middle ear (red dashed and blue solid curves 
respectively). These correspond to the vibration pattern 
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FIG. 1. BM response of the model compared to data sets in m taken at the (a-d) CF 6.6 kHz position and (e-h) at the 
position with CF 9.1 kHz. (a,e) BM displacement as a function of sound pressure level at different stimulus frequencies. (b,f) 
Frequency response functions at different SPLs, the lowest curve corresponds to 0 dB SPL. (c,g) Rate of growth indicating 
compression at different amplitudes as a function of frequency. (d,h) phase of the BM displacement. Thin black lines represent 
experimental data published in Rhode [21j and thick color lines represent numerical results. The bottom curve in (b-d) and 
(f-h) corresponds to 0 dB SPL, while consecutive curves (reading up) increase by 10 dB in SPL. 


of the most unstable mode of the cochlea. It can be seen 
that the magnitude of the vibration rapidly decreases to¬ 
wards the base of the cochlea when the middle ear is 
removed. However, the BM pattern is remarkably simi¬ 
lar to that with the intact middle ear for the rest of the 
cochlea. 

To explore the implications of local reflections, we also 
considered the case when the longitudinal coupling pa¬ 
rameter f q is increased by 20% at a single position, so 
that it forms a step function from base to apex, rather 
than varying randomly with length (data not shown). We 
might imagine that such an inhomogeneity will act as a 
single point of reflection. We found instead that distinct 
resonances were found (confirming the result in Ku et 
al. eg), and the cochlea patterns at the strongest res¬ 
onance were indistinguishable from the result in Figure 
[2jB,c). We therefore conclude that the theory of coherent 
reflections is not adequate to explain our results. 


D. Spontaneous emissions 

We can also use our unforced cochlear model to quan¬ 
tify the nature of the spontaneous oscillations that ex¬ 
ist beyond any point of instability. When an eigenvalue 
crosses the imaginary axis, as in the case of the rough 
unforced cochlea model described above, a Hopf bifurca¬ 
tion (of the entire cochlea) occurs; as a result, we expect 
a spontaneous oscillation to develop at approximately the 
frequency of the eigenvalue. Figure[3ja) shows the results 
of a combination of numerical continuation and simula¬ 
tion techniques to track both stable and unstable periodic 
orbits, and also to reveal more complicated dynamical be¬ 


havior. The resting position of the cochlea is represented 
by a blue line, the solid portion representing where this is 
stable and the dashed portion where it is unstable. The 
red curves show the results of tracking the limit cycle mo¬ 
tions born at a sequence of Hopf bifurcations that occur 
as the surface roughness a increases, shows a graph of 
the cochlear response as er varies. The vertical axis shows 
the stapes velocity v ow = x ow at the time instance when 
the stapes displacement is zero. This creates a sequence 
of values v ow which can be represented as a so-called 
Poincare map v ow K > P{v ow ). Superimposed on the plot 
are the results of direct simulations, after any transient 
motion has decayed. These are represented by black and 
green dots. The black dots were obtained by making a 
slow forward sweep in er and the green dots by making a 
backward sweep. Corresponding spectrograms are shown 
in Figure [3}'1>.<:). 

Note from Figure [^a) that as er increases from zero, 
the unforced cochlea goes unstable at around er = 0.05, 
via a supercritical Hopf bifurcation. Thus, a stable peri¬ 
odic orbit is born as er increases through this value. This 
represents a single frequency SOAE at «13.072 kHz, seen 
clearly in the spectrograms. Further supercritical Hopf 
bifurcations occur from the (now unstable) steady state 
as er increases further, although these branches of limit 
cycles all appear to be unstable. The single stable peri¬ 
odic solution can be inferred from where the where black 
dots of the simulation data overlie the red line from nu¬ 
merical continuation. The amplitude of this limit cycle 
grows steadily with er. 

Suddenly, for er « 0.13738, the (single frequency) pe¬ 
riodic orbit loses stability at a = 0.13738, through a su¬ 
percritical secondary Hopf (or Neimark-Sacker) bifurca- 
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Distance from stapes 


FIG. 2. (a) Characteristic roots of the cochlear model. Green 
circles correspond to the smooth cochlea, red crosses to the 
rough cochlea with er = 0.06, and blue pluses represent a 
cochlea of the same roughness without the middle ear. The 
horizontal axis is the frequency of the characteristic root and 
the vertical axis is the relative damping, which is the log of 
the amplitude ratio of two consecutive periods of free vibra¬ 
tion with the given frequency. (b,c) Magnitude and phase 
of the eigenvector corresponding to the most unstable eigen¬ 
value, for the rough cochlea model; red dashed and solid blue 
lines correspond to models with and without the middle ear 
respectively. 


tion. This creates quasi-periodic motion represented by 
the wide region of black dots in Figure |3ja). A second 
frequency 1.93 kHz is clearly visible in the spectrograms. 
At er ss 0.168 a second quasi-periodic orbit appears, with 
frequencies 3.16 and 13.06 kHz, which coexists with the 
first for larger values of a. The two different attractors 
are identified by different patterns of black and green dots 
in Figure [3ja) for 0.168 < a < 0.186, and different pat¬ 
terns in the spectrograms in Figure [3^b,c) . At a « 0.186 
the first quasi-periodic orbit (black dots) breaks up to 
form a chaotic attractor; this is indicated by a sudden 
broadening of the spectrogram in Figure j3jc). 

Figure |4]illustrates the qualitatively different classes of 
emissions predicted by the model in more detail, show¬ 
ing both the orbits (as Poincare maps) and their spectra. 
The first quasi-periodic orbit, just after onset, is shown 
in Figure [dja-b) and the two coexisting quasi-periodic 
orbits by green and black dots in Figure [IJc-d); note 



FIG. 3. (a) Amplitude of spontaneous emissions calculated 
from the model, represented by oval window fluid velocity as 
a function of surface roughness. Red (periodic orbits) and 
blue (fixed point) curves denote the results of numerical con¬ 
tinuation of the model. Black and green dots represent simu¬ 
lation results when increasing and decreasing cr, respectively. 
Points indicated by arrows represent responses plotted in de¬ 
tail in Figure [4] (b,c) Spectrograms of the cochlea emission 
with decreasing a and increasing er respectively. 


how the differing second frequencies result in different 
spacings between the peaks in the spectra. The chaotic 
attractor, shown with black dots in Figure |4je-f ) is of 
broadband character. Finally, the chaotic attractor dis¬ 
appears at er = 0.1892, where only a single quasi-periodic 
orbit remains, shown in Figure |4](g-h). 

The scenario described above can change significantly 
when applying different roughness patterns to the feed¬ 
back coefficient f q . However, simulations with different 
realizations (not shown) suggest that the supercritical na¬ 
ture of the bifurcation remains, and that supercriticality 
is the property of the chosen open probability function 
Pq . Changing Pq can result in a subcritical Hopf bifur- 
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FIG. 4. Poincare maps (a,c,e,g) and corresponding frequency spectra (b,d,f,lr) of solutions of the cochlear model, for values of 
cochlear roughness a indicated in Figure [3ja). See text for details. 


cation, as reported in m- 

The result that at a single parameter value multiple 
stable vibration patterns can coexist has implications 
for explaining other features observed experimentally. A 
large enough perturbation (e.g., a click) can push the 
ear into exhibiting different SOAE spectra, making one 
SOAE frequency appear and another disappear. Such 
changes are found experimentally [5j, as are hallmarks 
of chaotic oscillations [12], in qualitative agreement with 
the dynamics we predict. 


III. CONCLUSIONS 

In this paper we have introduced a nonlinear 
transmission-line model of the mammalian hearing or¬ 
gan. We believe this to be the first mathematical model 
that is both able to capture the turning curves of the 
Basilar membrane across a range of different frequen¬ 
cies and come up with credible, testable explanations for 
observed finite-amplitude otoacoustic emissions. Specif¬ 
ically, the model provides a parsimonious synthesis of 
many of the features that have been proved important 
in previous modeling and experimental studies over sev¬ 
eral decades: fluid-structure interaction leading to waves 
that travel a frequency dependent distance to the charac¬ 
teristic frequency position along the Basilar membrane; 
longitudinal variation in organ of Corti properties; local 
nonlinearity due to the combined effects of outer hair cell 


somatic motility and hair bundle adaptation; and spatial 
feed-forward effects due to longitudinal coupling within 
the OOC. 

The results in Figure [T] show how the combination of 
longitudinal coupling together with local nonlinearity is 
able to reproduce both qualitatively and quantitatively 
BM vibration data and neural tuning in different ani¬ 
mals. We are not aware that models with either pure 
longitudinal coupling or pure local nonlinearity are able 
to reproduce such features. 

Moreover, by applying a range of techniques from dy¬ 
namical systems theory, we have been able to calculate 
both instability thresholds and the waveforms of post¬ 
instability spontaneous otoacoustic emissions. We have 
shown how the roughness of the cochlea can be responsi¬ 
ble for producing these resonant instabilities. Moreover, 
by removing the middle ear and allowing no reflections 
at the oval window, we found that the resonance to only 
slightly diminish, but for the overall pattern of BM mo¬ 
tion to be maintained. This shows that SOAEs can be 
produced by instabilities that do not require reflections 
from the middle ear. Reflections, however, might occur 
elsewhere. 

Our model also how spontaneous OAEs arise robustly 
from supercritical Hopf bifurcations of the entire cochlea. 
Thus we do not need to make artificial “Hopf ear” hy¬ 
potheses about local nonlinearities at a particular CF 
being responsible. Furthermore, with significant rough¬ 
ness, the model exhibits a range of complex emissions be¬ 
yond pure tones including chaotic signals and bi-stability 



















































between different kinds of SOAE. These properties can where k = m ow /m, 8 = n/k ow and G e = G/m ow . 
explain reported recordings of emissions that appear to 
spontaneously switch between different signal patterns. 


IV. MATERIALS 

A. Equivalent formulation of the model 


B. Numerical methods 


In order to be able to discretize our model into a set of 
ordinary differential equations, we first combine ([T]) and 
© to eliminate the pressure p 


d 2 


d 2 




^bm £ *£bm — E ( ^C^0*tbm(£; t) AcJo^bm^? f) 


fqU'O 


J w(h/L)q (£ — h, t) d/i^ . (5) 


Equation ^ stays unmodified. The governing equations 
of the middle ear are rewritten in state-space form using 
X\ = 3Tm(0 ,t), x 2 = x hm (0,t), yi = x ow and y 3 = i ow - 
(m/m ow )xb m (0,t) as follows 


We discretize our equations along the length of the 
cochlea using finite differences. We use a non-uniform 
mesh that has interval length between the mesh points 
proportional to L(f). The first space derivative is a back¬ 
ward looking finite difference, and the second derivative 
is obtained using central differencing. The integral rep¬ 
resenting the feed-forward is approximated by the trape¬ 
zoid rule. Discretizing the governing equations with this 
scheme yields a set of ordinary differential equations that 
are solved using the MATLAB routine ode 113, specify¬ 
ing a constant mass matrix that arises from the second 
spatial derivatives of the right hand side of <§• The rel¬ 
ative error tolerance was set to 10 -10 , and the absolute 
tolerance to 10 -11 . 


Vi = 2/3 + k 1 2 3 4 5 6 7 8 9 10 ®2, 

V3 = — Cow^'ow !J 3 A K UJfo m Xi A G e [) e , 

±1 = x 2 , 


x 2 — 8 



8 rj (2Cbm^ J bm^'2 A 

(^Cbm^bm 2^ ow W ow ) X 2 UJfo m X 1 

K (G e p e 2(^ 0w CU 0 w 2/3 A>wl/l ) ) 


The numerical continuation to obtain both the pure- 
tone response in Figure [3] and the periodic orbits branch¬ 
ing from the Hopf bifurcation points are calculated from 
a periodic boundary value problem in time. The tem¬ 
poral discretization is performed using orthogonal collo¬ 
cation |2 with 4th degree interpolating polynomials on 
12 equidistantly spaced intervals. We used pseudo ar- 
clength continuation [1] to detect bifurcation points, and 
grow branches of solutions. 
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